Sinteni blokovi

Tokom evolucije, dolazilo je do razdvajanja vrsta. Kako su se vrste udaljavale dolazilo je do pojave izmena u DNK sekvenci. Neki regioni su mutirali vise, neki manje. Regioni DNK koji su kodirali najbitnije proteine za opstanak jedinke nisu trpeli velike izmene (u protivnom jedinka ne bi prezivela), ali su se desavale izmene u redosledu navedenih regiona.

Primer: Neka je redosled određenih regiona kod jedne vrste: $A \rightarrow B \rightarrow C \rightarrow D$

Kod druge vrste, koja se od zajedničkog pretka sa prvom vrstom redosled istih regiona je: $D \rightarrow C \rightarrow A \rightarrow B$

Između regiona postoje sekvence proizvoljne

Sekvence u regionima mogu biti inverzne (sekvenca se cita u obrnutom redosledu u odnosu na pretka): $(-D) \rightarrow C \rightarrow (-A) \rightarrow B$

Preciznije, neka su regioni A = ATAG, B = CCGT, C = GCTG, D = TTTA

Jedna vrsta: $A \rightarrow B \rightarrow C \rightarrow D$ = ATAG...CCGT...GCTG...TTTA

Druga vrsta: $(-D) \rightarrow C \rightarrow (-A) \rightarrow B$ = ATTT...GCTG...ATAG...CCGT

Regioni koji se posmatraju u ovom kontekstu se nazivaju sinteni blokovi (eng Synteny Blocks). Razlike u redosledu i pravcu citanja sinteni blokova nastaju usled savijanja i presecanja DNK. Pojedinačne promene koje dovode do inverzije redosleda sinteni blokova nazivaju se inverzije.

Primer: $A \rightarrow | B \rightarrow C \rightarrow D |$ Posle inverzije u označenom regionu: $A \rightarrow (-D) \rightarrow (-C) \rightarrow (-B)$

Broj inverzija između dve vrste često je srazmeran periodu koji je protekao od odvajanja vrsta od zajedničkog pretka do sadašnjosti.

Jedan od pristupa za računanje broja inverzija je, naravno, pristup grubom silom. Kako je ovaj pristup veoma neefikasan, može se probati sa heurističkim, znatno bržim pristupom. Jedan takav, heuristički, pristup je Sorting by Reversals. Sinteni blokovi jedne vrste označe se redom brojevima od 1 do $n$ (gde je $n$ ukupan broj blokova) dok se odgovarajući blokovi druge vrste označavaju oznakama kojim su označeni u prvoj vrsti.

Primer: Prva vrsta: $A \rightarrow D \rightarrow C \rightarrow B$ = $[1, 2, 3, 4,]$

Druga vrsta: $D \rightarrow (-A) \rightarrow C \rightarrow (-B)$ = $[2, -1, 3, -4,]$

Polazeći od drugog, neuređenog niza, primenom inverzija se niz dovodi do uređenog niza poziciju po poziciju.

Neka je početno stanje: $[2, -1, 3, -4,]$, kako je na prvoj poziciji pogrešan blok, potrebno je na prvu poziciju dovesti blok broj 1.

$[2^*, -1, 3, -4,]$

$[|2, -1|, 3, -4,]$

$[|1, -2|, 3, -4,]$

Ukupni broj operacija potrebnih za dovođenje nesortiranog niza do sortiranog predstavlja broj inverzija između dve vrste. Ovo je, kao što je napomenuto, heuristički pristup koji ne vraća uvek najmanji broj inverzija. Kako je očekivano da je evolucija ipak išla najkraćim putem, ima smisla tragati i za boljim rešenjem.

In [1]:
def apply_sorting_reversal(P, k):
    n = len(P)
    for i in range(k, n):
        if abs(P[i]) == k + 1:
            L = P[:k]
            Mid = P[k: i+1]
            Mid = [-x for x in Mid[::-1]]
            R = P[i + 1: n]
            
            return L + Mid + R
In [2]:
def greedy_sort(P):
    approx_reversal_distance = 0
    
    # Optional
    print(P)
    
    n = len(P)
    for k in range(n):
        if P[k] != (k + 1):
            P = apply_sorting_reversal(P, k)
            approx_reversal_distance += 1
            print(P)
            
            if P[k] == -(k + 1):
                P[k] = -P[k]
                approx_reversal_distance += 1
                print(P)
                
    return approx_reversal_distance
In [3]:
P = [+1, -7, +6, -10, +9, -8, +2, -11, -3, +5, +4]
greedy_sort(P)

[1, -7, 6, -10, 9, -8, 2, -11, -3, 5, 4]
[1, -2, 8, -9, 10, -6, 7, -11, -3, 5, 4]
[1, 2, 8, -9, 10, -6, 7, -11, -3, 5, 4]
[1, 2, 3, 11, -7, 6, -10, 9, -8, 5, 4]
[1, 2, 3, -4, -5, 8, -9, 10, -6, 7, -11]
[1, 2, 3, 4, -5, 8, -9, 10, -6, 7, -11]
[1, 2, 3, 4, 5, 8, -9, 10, -6, 7, -11]
[1, 2, 3, 4, 5, 6, -10, 9, -8, 7, -11]
[1, 2, 3, 4, 5, 6, -7, 8, -9, 10, -11]
[1, 2, 3, 4, 5, 6, 7, 8, -9, 10, -11]
[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, -11]
[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11]
Out[3]:
11

Drugi pristup dolazi iz sveta grafova. Niz sinteni blokova predstavi se kao cikličan, umesto do sada linearne pretpostavke. Jedan sinteni blok predstavlja se usmerenom granom grafa između dva čvora, simbolizujući pravac čitanja bloka. Takođe, kako je moguće da se blokovi nalaze u različitim hromozomima, formiraju se posebni ciklusi za blokove svakog hromozoma.

Grane koje predstavljaju blokove oznacene su crnom bojom dok su grane između blokova (koje predstavljaju sekvence između sinteni blokova) označavaju obojenim granama.

Jedna inverzija na hromozomu rezultovaće dvostrukim ukrštanjem grana (2 break) grafa, na način prikazan na slici: Preseci koji predstavljaju inverzije mogu se desiti samo nad obojenim granama.

Da bi se izračunao broj inverzija između dve vrste (grupe hromozoma), formira se zajednički, Breakpoint graf, preklapanjem grana obe grupe. Zatim, računa se broj potrebnih inverzija kako bi se dobijeni graf "razmrsio", tj. doveo do grafa koji sadrži samo proste cikluse.

Traženi broj inverzija dobija se uzastopnim pravljenjem preseka nad Breakpoint grafom, tačnije, presecanjem grana jedne boje (npr. samo crvene, gde su crvenim grana predstaljeni hromozomi jedne vrste) koje pripadaju ciklusima sa granama druge boje, sve dok graf ne postane prost (sve obojene grane se poklapaju između istih čvorova).

Naredni kodovi odgovaraju pseudokodovima iz udžbenika i primeri su postavljeni da odgovaraju primerima sa slika.

In [4]:
def chromosome_to_cycle(chromosome):
    n = len(chromosome)
    nodes = [0] * (n*2)
    for j in range(n):
        i = chromosome[j]
        
        if i > 0:
            nodes[2 * j] = (2 * i) - 1
            nodes[2 * j + 1] = (2 * i)
        else:
            nodes[2 * j] = - (2 * i)
            nodes[2 * j + 1] = - (2 * i) - 1
            
    return nodes
In [5]:
chromosome_to_cycle([1, 2, 3, 4, 5, 6])
Out[5]:
[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12]
In [6]:
def cycle_to_chromosome(nodes):
    n = len(nodes)
    chromosome = [0] * (n//2)
    for j in range(0, n, 2):
        if nodes[j] < nodes[j + 1]:
            chromosome[j // 2] = nodes[j + 1] // 2
        else:
            chromosome[j // 2] = - nodes[j] // 2
            
    return chromosome
In [7]:
cycle_to_chromosome([1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12])
Out[7]:
[1, 2, 3, 4, 5, 6]
In [8]:
def colored_edges(P):
    edges = []
    for chromosome in P:
        nodes = chromosome_to_cycle(chromosome)
        n = len(nodes)
        for j in range(1, n - 1, 2):
            edges.append((nodes[j], nodes[j + 1]))
        edges.append((nodes[n - 1], nodes[0]))
            
    return edges
In [9]:
P = [[+1, -2, -3, +4]]

nodes = chromosome_to_cycle(P[0])
print(nodes)
chromosome = cycle_to_chromosome(nodes)
print(chromosome)
P_colored_edges = colored_edges(P)
print(P_colored_edges)

[1, 2, 4, 3, 6, 5, 7, 8]
[1, -2, -3, 4]
[(2, 4), (3, 6), (5, 7), (8, 1)]
In [10]:
def black_edges(P):
    edges = []
    for chromosome in P:
        nodes = chromosome_to_cycle(chromosome)
        n = len(nodes)
        for j in range(0, n - 1, 2):
            edges.append((nodes[j], nodes[j + 1]))
            
    return edges
In [11]:
print(black_edges(P))

[(1, 2), (4, 3), (6, 5), (7, 8)]

Klasa koja predstavlja graf

In [12]:
class GenomeGraph():
    def __init__(self, edges):
        G = {}
        
        for edge in edges:
            (v, w) = edge
            
            if v not in G:
                G[v] = [w]
            else:
                G[v].append(w)
                
            if w not in G:
                G[w] = [v]
            else:
                G[w].append(v)
                
        self.G = G
                
    def get_cycles(self):
        unvisited = set(self.G.keys())
        
        cycles = []
        
        while len(unvisited) > 0:
            start = min(unvisited)
            v = start
            unvisited.remove(start)
            
            cycle = []
            
            while True:
                next_w = None
                for w in self.G[v]:
                    if w in unvisited:
                        next_w = w
                        break
                            
                if next_w == None:
                    cycle.append((v,start))
                    break
                
                cycle.append((v, next_w))
               
                v = w
                unvisited.remove(w)
            
            cycles.append(cycle)
                    
        return cycles
    
    def remove_edge(self, edge):
        (v, w) = edge
        self.G[v].remove(w)
        self.G[w].remove(v)
        
    def add_edge(self, edge):
        (v, w) = edge
        self.G[v].append(w)
        self.G[w].append(v)
In [13]:
P = [[1,-2,-3,4], [5,6]]
genome_graph = GenomeGraph(black_edges(P) + colored_edges(P))
genome_graph.get_cycles()
Out[13]:
[[(1, 2), (2, 4), (4, 3), (3, 6), (6, 5), (5, 7), (7, 8), (8, 1)],
 [(9, 10), (10, 11), (11, 12), (12, 9)]]
In [14]:
def nodes_on_path(path):
    nodes = []
    
    for (v, w) in path:
        nodes.append(v)
    
    return nodes

Pretvaranje hromozoma u graf

In [15]:
print (P)
genome_graph = GenomeGraph(black_edges(P) + colored_edges(P))
cycles = genome_graph.get_cycles()
print(cycles)

for cycle in cycles:
    print(nodes_on_path(cycle))

[[1, -2, -3, 4], [5, 6]]
[[(1, 2), (2, 4), (4, 3), (3, 6), (6, 5), (5, 7), (7, 8), (8, 1)], [(9, 10), (10, 11), (11, 12), (12, 9)]]
[1, 2, 4, 3, 6, 5, 7, 8]
[9, 10, 11, 12]
In [16]:
def graph_to_genome(genome_graph):
    P = []
    cycles = genome_graph.get_cycles()
    
    for cycle in cycles:
        nodes = nodes_on_path(cycle)
        chromosome = cycle_to_chromosome(nodes)
        P.append(chromosome)
        
    return P
In [17]:
genome_graph = GenomeGraph(black_edges(P) + colored_edges(P))
graph_to_genome(genome_graph)
Out[17]:
[[1, -2, -3, 4], [5, 6]]

2-break nad grafom:

In [18]:
def two_break_on_genome_graph(genome_graph, i, i_p, j, j_p):
    genome_graph.remove_edge((i, i_p))
    genome_graph.remove_edge((j, j_p))
    
    genome_graph.add_edge((i, j))
    genome_graph.add_edge((i_p, j_p))
    
    return genome_graph
In [19]:
P = [[1, -2, -3, 4]]

genome_graph = GenomeGraph(black_edges(P) + colored_edges(P))
print(genome_graph.get_cycles())

genome_graph = two_break_on_genome_graph(genome_graph, 1, 8, 3, 6)
print(genome_graph.get_cycles())

[[(1, 2), (2, 4), (4, 3), (3, 6), (6, 5), (5, 7), (7, 8), (8, 1)]]
[[(1, 2), (2, 4), (4, 3), (3, 1)], [(5, 6), (6, 8), (8, 7), (7, 5)]]
In [20]:
def two_break_on_genome(P, i, i_p, j, j_p):
    genome_graph = GenomeGraph(black_edges(P) + colored_edges(P))
    genome_graph = two_break_on_genome_graph(genome_graph, i, i_p, j, j_p)
    P = graph_to_genome(genome_graph)
    return P
In [21]:
P = [[1, -2, -3, 4]]
print(P)
P = two_break_on_genome(P, 1, 8, 3, 6)
print(P)

[[1, -2, -3, 4]]
[[1, -2], [3, -4]]
In [22]:
P = [[1, 2, 3, 4]]
print(P)
genome_graph = GenomeGraph(black_edges(P) + colored_edges(P))
genome_graph = two_break_on_genome_graph(genome_graph, 2, 3, 4, 5)
genome_graph = two_break_on_genome_graph(genome_graph, 5, 3, 7, 6)
P = graph_to_genome(genome_graph)
print(P)

[[1, 2, 3, 4]]
[[1, -2, -3, 4]]
In [23]:
import random

def shorthest_rearangement_scenario(P, Q):
    print(P)
    red_edges = colored_edges(P)
    blue_edges = colored_edges(Q)
    
    breakpoint_graph = GenomeGraph(red_edges + blue_edges)
    
    while True:
        cycles = breakpoint_graph.get_cycles()
        selected_cycle = None
    
        for cycle in cycles:
            n = len(cycle)
            if n > 2:
                selected_cycle = cycle
                break
         
        if selected_cycle == None:
            break
        
        k = random.randrange(0, n)
        
        edge = selected_cycle[k]
        (v, w) = edge
        
        if (v, w) not in blue_edges:
            (w, v) = (v, w)
            
        if (v, w) not in blue_edges:
            k = (k + 1) % n
            
        (j, i_p) = (v, w)
        left_red_edge = selected_cycle[(k - 1) % n]
        right_red_edge = selected_cycle[(k + 1) % n]
        
        (i, j) = left_red_edge
        (i_p, j_p) = right_red_edge
        
        breakpoint_graph.remove_edge(left_red_edge)
        breakpoint_graph.remove_edge(right_red_edge)
        
        breakpoint_graph.add_edge((j, i_p))
        breakpoint_graph.add_edge((i, j_p))
        
        P = two_break_on_genome(P, j, i, i_p, j_p)
        print(P)
In [24]:
P = [[1, -2, -3, 4]]
Q = [[1, 2, 3, 4]]

shorthest_rearangement_scenario(P, Q)

[[1, -2, -3, 4]]
[[1, 2, -3, 4]]
[[1, 2, 3, 4]]
In [25]:
P = [[-6, 1, -2, 3, 4, 5]]
Q = [[1, 2, 3, 4, 5, 6]]

shorthest_rearangement_scenario(P, Q)

[[-6, 1, -2, 3, 4, 5]]
[[1, -2, 3, 4, 5, 6]]
[[1, 2, 3, 4, 5, 6]]